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Abstract 

We examine the effect of local matter on the chaotic behavior of a relativistic 
test particle in non-vacuum static axisymmetric spacetimes. We find that the sign 
of the sectional curvature in the geodesic deviation equation defined by the Riemann 
curvature does not always become a good tool to judge the occurrence of chaos in the 
non-vacuum case. However, we show that the locally unstable region ( LU region ) 
defined by the Weyl curvature can provide information about chaos even in non- vacuum 
spacetime as well as in vacuum spacetime. Since the Weyl tensor affects only the 
shear part of the geodesic congruence, it works effectively to stretch some directions of 
geodesic congruence, which helps to cause the chaotic behavior of geodesies. Actually, 
the orbit moving around an unstable periodic orbit (UPO) becomes strongly chaotic 
if it passes through an LU region, which means that the LU region can be used as a 
good tool to know in which situation the chaos by homoclinic mixing occurs around a 
UPO. 
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1 Introduction 

General Relativity (GR) is the most plausible theory of gravitation to explain relativistic 
astrophysical phenomena in our universe. It has succeeded in explaining many observational 
results in cosmology and astronomy by using idealized models such as Friedmann Robertson- 
Walker (FRW) spacetime or Schwarzschild spacetime. 

However these kinds of idealized models cannot cover all of the astrophysical phenomena 
which might actually happen in our universe. Chaos may be one of these phenomena because 
these idealizations usually force a strong symmetry on spacetime, which extinguishes the 
possibility of complicated behavior such as chaos. In a Hamiltonian system chaos can occur 
only in the non-integrable system where the number of integrals is less than the dimension 
of configuration space 0. 

For example, a Bianchi-IX universe, one of the spacetimes which is generalized from 
a closed FRW universe by dropping isotropy, is known to show stochastic behavior ||. 
However, there has been much discussion about the definition of chaos independent of the 
choice of time coordinate |J. 

As for test particle motion in GR, the Kerr-Newman spacetime, which includes the 
Schwarzschild, Kerr and Reissner-Nordstrom spacetimes as special cases, is known to be 
integrable. This can be explained by the fact that it is a Petrov type D spacetime and 
in addition to the energy E, angular momentum L and super-hamiltonian H = ^g^p^Pu, 
there is a fourth conservative quantity described by the Killing tensor J|. Thus bound orbits 
in this spacetime are strictly constrained to a torus in the phase space, and none of them 
behaves chaotically. However, in realistic situations the spacetime around compact objects 
such as a black hole or a neutron star is not expected to have such a high degree of symmetry, 
because of the distortion of the gravitational source itself or of the effect of other gravitational 
sources such as gravitational waves, magnetic fields or other astronomical objects. There is 
a possibility that the deviations of these spacetimes from Petrov D extinguish the Killing 
tensor and cause strong chaos in the test particle motion around those compact objects |5fl. 

In fact more generic spacetimes are usually non-integrable and some of them exhibit 
chaotic behavior in test particle motions |5|-|12||. In GR, an unstable periodic orbit (UPO) 



exists in axisymmetric and static spacetime in contrast with Newtonian mechanics. This is 
one of the main factors which causes chaos for geodesies in the spacetime. For example, chaos 
is caused by periodic perturbations around a UPO in an integrable static and axisymmetric 
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spacetime, such as the Schwarzschild spacetime fTHH [jTI]l . Although Schwarzschild spacetime is 
static, spherically symmetric and integrable, the perturbations extinguish the time symmetry 
or axisymmetry. This type of chaos can be explained by the homoclinic tangle around 



a UPO [13]. Such chaos can also occur by non-periodic perturbations like adding a small 
gravitational source to the spacetime, because it brings about the asymmetry around a UPO 
H [[J . On the other hand, it is not trivial to determine whether or not chaos occurs if there 
is no perturbation around a UPO. It seems necessary to find tools to judge where and with 
what strength chaos occurs in these spacetimes, in order to make a quantitative analysis of 
chaos. 

In the previous paper (paper I ||), we paid attention to the fact that in GR, free particle 
motion in a spacetime is described by a geodesic, so that local instability may be determined 
from the sign of the curvature tensor. We showed that in axisymmetric static vacuum 
spacetimes, geodesies become strongly chaotic after passing through an LU region defined 
by the distribution of eigenvalues of the Weyl tensor. Recently, Vieira and Letelier showed, 
with their careful analysis around UPO, that chaos predicted by the LU region also originates 
from a homoclinic tangle, so that the local instability criterion is not so strict as to predict 



the chaotic property for the orbit passing through an LU region fl4|| . This result suggests 
that in addition to the passage through an LU region, it seems necessary for the orbit to pass 
by the UPO to cause strong chaos. However, from a physical viewpoint, the origin of the 
chaos around a UPO is still unclear without information about the LU region. For example, 
why the torus suddenly begins to break around UPO in Fig. 5 of [Q] cannot be explained just 
by the homoclinic tangle. Since the LU region can be utilized as an effective tool for judging 
the occurrence of chaos even in such a case, it is still important to analyze the correlation 
between the LU region and chaos, and to give its physical meaning. 

In the previous paper, we examined the eigenvalues of the Weyl tensor, because in vacuum 
spacetime the Weyl tensor coincides with the Riemann tensor which is used to determine 
the local instability in the geodesic deviation equation. In the non-vacuum case, however, 
the situation is not simple, because the Riemann tensor is no longer the same as the Weyl 
tensor. In this paper, to see the effect of matter, we first analyze the deviation of geodesies by 
dividing the Riemann tensor into the Weyl tensor, the Ricci tensor and the scalar curvature 
in section 2. This procedure makes it clear how the eigenvalues of the Riemann tensor relate 
to those of the Weyl tensor. In section 3, we will numerically show that the LU region 
defined by the Weyl tensor still remains a good tool to determine chaos even in spacetimes 
with a matter field. We will also show that even if the local matter effect is strong enough 
to make the sectional curvature K(u,n) negative everywhere in bound regions, chaos can 
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still occur around an LU region. This result supports the claim of Vieira and Letelier on the 
weakness of the local instability determined by the sign of JC(u,n). In section 4, in order 
to reveal the physical origin of such chaos, we will analyze this type of chaos from another 
viewpoint, the shear effect of the geodesic congruence, and show the role of the Weyl tensor 
on the eigenvalues of the shear matrix. We will show that the LU region has the character of 
possessing two independent stretch directions for geodesic congruence, which might explain 
the contribution of the LU region to the chaos of geodesies passing by a UPO. Finally, we 
will give our conclusions and some remarks in section 5. 

2 Eigenvalues of the Riemann tensor 

2.1 Derivation of the eigenvalues of the Riemann tensor 

In this paper, we use the same notation as in paper I |§ . In GR, the motion of a free particle 
is described by a geodesic and its deviation n M is given by the equation 

— = -Rft^nW, (2.1) 

where is the 4-velocity and is orthogonal to u M , i.e., u^u^ = — 1, and u^n^ = 0. 
Equation (|2.1|) is also reduced by the sectional curvature JC(u, n)(= —R(u,n,u,n)/\\n\\ 2 ) of 
the plane spanned by u and n to, 



2 



^H=/C( tt) n)H + - 



Dn 

n x 

Dt 



2 



(2.2) 



where we use the notation ||V|| = (| V^V^) 1 ' 2 for a 4-vector V^. From Q2.2|) , if fC(u,n) is 
positive, a geodesic becomes locally unstable since the second term in the right hand side of 
( |2.2|) is always positive. As described in paper I, the sectional curvature takes critical values 
in the direction of the eigenvectors of the Riemann tensor. Here we will show how these 
eigenvalues are expressed in a static axisymmetric spacetime. 

In general, the Riemann tensor R^ upa can be decomposed by using the Weyl tensor C^ pa , 
the Ricci tensor R pv and the Ricci scalar R into the form 

Rpupo = C^upu + Qp\pR<j\v 9u[pR<j]p T^Rff n[p9(r]u ■ (2-3) 

In bivector formalism, the Riemann and Weyl tensors are regarded as 6 x 6 matrices which 
are decomposed as follows: 
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where £ , £, H, H, T and T are 3x3 matrices |T(|. In the vacuum case the parts composed 
of the Ricci tensor and Ricci scalar vanish from Einstein's equation, and the Riemann tensor 
coincides with the Weyl tensor. However, in the non-vacuum case, we do not have such an 
advantage. We have to deal with both 1Z and C, separately. If the spacetime is static, the 
matrices 7Z and C are expressed as follows: 

In this case, we can always diagonalize the matrices TZ and C by using an appropriate 
orthonormal real tetrad basis {e( ), e^, e( 2 ), e( 3 )} since the matrices £, £ and T in 71 and C 
are all symmetric. We denote the six eigenvalues of 7Z as follows 

K oi = -^(e(o),e w ,e { o),e (i) ) (i < j, i, j = 1 ~ 3) 



i?(e (i) , e {fc ), e (i) , e (fc) ) = (j < fc, j, fc = 1 ~ 3). (2.6) 



Then, the sectional curvature fC(u,n) is expressed by k£ v as follows. 

K.(u,n) = —R(u,n,u,n) = —u^n^u^n^R(e( a ),e(i3),e(j),e(€)) 

= ^S^(uM n M -uM n M)\* (2.7) 



where 5^ is 1 if v is equal to zero, and otherwise it is —1. From $L% , a positive value of 
contributes to the local instability ( K.(u,n) > 0), while a positive value of contributes 
to the local stability ( K,(u,n) < 0). 

For static and axisymmetric spacetimes, the metric is written as 

where t and <fi are the coordinates related to two Killing vectors d/dt and d/dcj) and both U 
and /c are functions depending only on p and z |. In these coordinates, e(o) and e(3) coincide 
with the normalized Killing vectors eg and eg, respectively. Then and K,f 2 always become 
7Z 3 3 and 7Z 6 6 , respectively and the other eigenvalues are 

1 



< " 2 
K ° 2 ~ 2 



2 \2 



(^ 2 2 + 7e 1 1 )- v /(^ 2 2 -^ 1 1 ) 2 + 4(^ 2 1 ) 2 



(2.9) 



where the matrix 7Z A B is the Riemann tensor R in bivector formalism and the components 
A and B of 7^ run 1 ~ 6 which denote the tetrad components (ip), (tz), (i(j>), (z(j>), (4>p), 
* We use units of G = c = 1, but we explicitly write G or c when it may help our discussion 
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and (pz), respectively. and be derived by changing the indices 1 and 2 in 

to 4 and 5, respectively. We can also get the eigenvalues of C using the same tetrad basis as 
follows, 



Of 

c 



jk 



-C(e( ), e(i), e (0 ), e^) 
C( e (i)j e (fc) ? e (i)> e (fc)) : 



(j < fe, j, k 



3). 



(2.10) 



Here we have 



For the metric (J2.8|) , these eigenvalues, k° are derived by changing 



the Riemann tensor 7Z in 



K 
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K. 
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into the Weyl tensor C as 



Cl + 



d) 2 + ±{c\y 



(2.11) 



and = k\ 2 = C\. In [§], we called a region where k,q 3 < 0, Kq 2 > and Kq X > a locally 
unstable (LU) region, since > and Kq 2 > means that the sectional curvatures spanned 
by e(o) and one of the two eigenvectors orthogonal to two Killing directions, /C(e(o), e(i)) and 
/C(e(o), e(2)) are positive. In fact, we showed in || that this region strongly affects the chaotic 
behavior of free particle motion. 

In the non-vacuum case, however, the eigenvalues of the Riemann tensor are not the 
same as those of the Weyl tensor and no longer degenerate. We cannot divide the spacetime 
by the eigenvalues of Riemann tensor as simply as in the vacuum case, because we have six 
independent eigenvalues. 

In order to judge the chaos from curvature information, we have to know the relationship 
among these eigenvalues and clarify how the effect of local matter changes the value of k c ^ v 
from that of k^ u . By using the relationship between the Riemann tensor and the Weyl tensor 
( [2.3D and Einstein's equation = 87rT M1/ , we get the relationship between k c v and n^ u (see 
Appendix A for the detail). 



v 0» 



K 
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jk 
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T 
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{po -Pi + 2(pj + p k )} 



-jk + y{5(po 



Pi) - 2(Pj+Pk)}, 



(2.12) 



where p = T(e {0 ), e {0 )) and p» = T(e^, e 



(0/ 



From (|2.7| ), both positive values of and negative values of nf k contribute to make 
JC(u, n) positive. In general we cannot give any concrete inequality between k° and by 
specifying energy conditions. 
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However, we can find that the electromagnetic field always weakens the local instability 
determined by the Weyl curvature. That is, in this case T vanishes, so that ( |A5| ) in Appendix 
A becomes 



"Oi 
.R 



«w -4tt{po -Pi} 



K? k = K c jk + 4n{p - Pl } (i,j,k = 1~3) (2.13) 
and po ~ Pi is always positive, because of the dominant energy condition [ITH]. Hence the 



eigenvalues are always smaller than k^, while nf- is always bigger than 

If the matter effect is not strong enough to change the sign of Kq 2 or Kq 3 in the LU region, 
the region where both k^ 2 and Kq 3 are positive, which we call the RLU region, remains inside 
the LU region. Then we may use the RLU region for a criterion of chaos. However, if the 
RLU region does not exist inside the LU region, it seems to be difficult to explain chaos 
around the LU region by the positivity of K.(u, n). In fact in the vacuum case, the LU region 
was characterized as the region where /C(e(o), e(i)) > and JC(e(o), &{2)) > in ||. However 
it is nothing but the RLU region, if a matter field exists. 

In the next section, we will numerically examine the relationship between the chaotic be- 
havior of a geodesic and its passage through an LU or RLU region by using some analytically 
given spacetimes. 



3 Numerical analysis 

Here we examine the relationship between chaos and LU or RLU regions by numerical 
calculations for some spacetimes with matter fields. First, we reanalyze the Majumdar- 
Papapetrou(MP) spacetimes as examples of spacetimes where both LU and RLU regions 
appear. Then, we study the Ernst universe as an example of a spacetime where no RLU 
region exists inside the LU region. Finally, we treat the spacetime with a scalar field as an 
example of the spacetime with a matter field which does not satisfy any energy condition. 



3.1 Chaos around N-maximally charged black holes 

First, we examine the MP solution. This describes the system of several extreme Reissner- 
Nordstrom (RN) black holes that are balanced by the attractive gravitational and the re- 
pulsive electric forces. In this spacetime, N black holes can be located at random in a three 
dimensional hypersurface. The metric is 



ds 2 = -V- 2 dt 2 + V 2 (dx 2 + dy 2 + dz 2 ) (3.14) 
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N M, 



v = i + > ( 3 - 15 ) 



i=l r 



where Mj is the mass of the zth black hole and = [(x — Xi) 2 + (y — yi) 2 + (z — Zj) 2 ] 1 / 2 is the 
Euclidean coordinate distance from the location of the ith black hole, (xi,yi,Zi), to a point 
(x,y,z). If iV = 1, it is just the extreme RN solution. Since the RN solution belongs to 
Petrov type D, two eigenvalues are degenerate everywhere and neither LU nor RLU regions 
appear. This is consistent with the fact that stationary axisymmetric Petrov D spacetimes 
are integrable and chaos does not occur there. 



The case of N > 1 has been examined by several authors 0112]. They indicated that 
these spacetimes cause chaotic type motions characterized by fractal boundaries of initial 
conditions for both timelike and null geodesies on the meridian plane. 

Here we examine the case with the N black holes located on the z axis, since we are 
interested only in the axisymmetric case. In this case, the spacetime is described by using 
the two functions U and k in the metric (12.81) with 



U = -\n{l-U c ), 

k = 0, (3.16) 



where Uc corresponds to that of iV-Curzon spacetime [2C] and is 



N M- 

Efc = -<?£—■ (3-17) 



i=l 1 



First, we study the case N = 2 with non-zero angular momentum of a test particle, L. In 
this case, we can easily find an LU region as in 2-Curzon or 2-ZV spacetime ||, since each of 
two gravitational sources on the z axis attracts the geodesies toward itself (Fig.|l|(a)). Then 
if the energy of test particle is appropriately large, the bound orbit can be overlapped with 
LU region and strong chaos occurs beyond some critical energy (Fig.[TJ(b)). Here we can also 
find a RLU region inside an LU region, since the effect of the electric field produced by the 
electric charge of each black hole is not so strong as to extinguish it. 

Next, we consider the case of N — 3, with one of the black holes located at the origin and 
the others at the same distance on both sides of the z-axis. In this case, we can also see an 
RLU region inside an LU region. When the mass of the central black hole is much more than 
that of the others, the LU region does not overlap the bound region. We have not seen any 
chaotic behavior there, even if we increase the energy up to the value of the UPO. However, 
as the mass of the central black hole reduces to that of the others, the bound region are 
more and more overlapped with the LU region. Here we analyzed the case that the energy 
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of geodesic is a little higher than the energy of UPO, Eupo, where the classically admitted 
region does not become a bound region, but has the tiny throat connecting it with a black 
hole. In this case, we can check whether or not the chaos determined by an LU region always 
has the origin of homoclinic mixing as Vieira and Letelier claimed, since the orbit passing by 
the throat falls into a black hole. In this case when the orbit passes through the LU region, 
the torus in Poincare map is completely broken(Fig.^), although the orbit falls into black 
hole with a finite time interval. These results certainly suggest that passing through the LU 
region (and RLU region) contributes to the occurrence of chaos (though it must be defined 
with finite time interval) as for geodesic motion in a multi-black hole system, regardless of 
homoclinic tangle. 



3.2 Chaos around the Ernst universe 

Next, we reexamine the Ernst solution which includes a magnetic field. Axisymmetric sta- 
tionary spacetimes with magnetic field B Q along the z axis can be derived from the corre- 
sponding vacuum spacetimes In particular, when we start with the metric form 



for the static vacuum case, we find the function Ub and ks of the spacetime with magnetic 
field to be : 

U B = f/ln(l + ii? Ve- 2C/ ) 

k B = k + 2\n(l + ^Bye- 2U ), (3.18) 

where B is the magnitude of the magnetic field along the z axis. 

Several years ago, Karas and Vokrouhlicky found that both neutral and charged particles 
behave chaotically in one of the Ernst solutions, i.e, around a Schwarzschild black hole 
immersed in a magnetic field JF|. These particles show chaotic behavior more and more 
strongly as the energy increases and approaches the critical value, Eupo, beyond which the 
particles fall into the black hole. This suggests that the chaos occurs through the homoclinic 
mixing around a UPO. However the origin of such a homoclinic mixing is still unclear, 
at least in the respect that the metric has a reflection symmetry on the equatorial plane to 
retain a homoclinic orbit around the UPO. Here we examine whether it is possible to explain 
the chaos by passing through an LU region. We reanalyze the case which they examined, 
that is, the case that U and k correspond to the functions of Schwarzschild spacetime and 
Bq = 0.15 /M. In Fig.||] we plot the LU region and bound region of a neutral test particle 
with (E 2 , L) = (29.12 fi 2 , 25.0 fiM), where \i is the rest mass of a test particle. From this 
figure, we can see that the magnetic field B plays a role in causing an LU region around 
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a single black hole, since no LU region appears for the case of So — 0. Certainly, it is 
possible to explain the strong chaos shown in |F] by the passage through an LU region, 
since the orbit moves inside the bound region almost ergodically and passes through the 
LU region regardless of the initial conditions. On the other hand, we cannot see an RLU 
region anywhere inside an LU region, because the matter effect is strong enough to change 
the sign of Kq X or Kq 2 into negative through equation ( [2.13| ). This result makes it difficult to 
explain chaos by local instability determined by geodesic deviation equation, since the local 
instability in the p-z plane is determined by and Kq 2 , rather than «§i an d ^oV 

So here, in order to see the effect of the absence of an RLU region, we consider the 
eigenvalues of the following matrices, 

t% = C» ap u»u°. (3.19) 

The sectional curvature, JC(u, n) for any pair of (u, n) is described by the linear combination 
of the non-zero eigenvalues a-, (i — 1, . . . , 3) of matrix 71 as follows, 

3 

K(u,n) = £)a s (rjf. (3.20) 

i=l 

(See Appendix B.) In order to fix the 4-velocity u in ( |3.19 ), we use the parameter 0* defined 

as, 

= V* cos 0* 

w (2) = ^sin©„ (3.21) 
where v * is the meridian velocity which we defined in M as follows. 

v,{x,E,L)- 2| , Q/£UI |, , (■>.>>) 



where x = (p, z) and is the effective potential for the particle with the angular momentum 
L. The rest components, and are determined by E and L, respectively. 

We numerically examined the distribution of eigenvalues a\ in the bound region and 
found that all of the eigenvalues a- become negative everywhere inside the bound region for 
any pair of (f*,0*) with (E 2 , L) = (29.12 /i 2 , 25.0 pM). This means that chaos cannot 
be explained by the curvature term in Q2.1|) or Q2.2|) , since K(u,n) becomes negative for any 
pair of (u,n) from ( 3.20|) . This result seems inconsistent with the result that strong chaos 
occurs for the same value of E and L, which leads us to new analysis by a shear effect of 
the geodesic congruence on chaos rather than that of the geodesic deviation, as we will see 
in section 4. 
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3.3 The spacetime with scalar field 

Finally, we examine a spacetime with a scalar field, as an example of a spacetime with a 
matter field which satisfies neither the strong nor dominant energy condition. In this case, 
an RLU region need not be included in the LU region, since the matter field may increase 
the local instability. 

Here we make use of the solution in which iV-scalar charged ZV-type singularities are 



located on the z axis [|2T|. We examine the case of two scalar charged singularities balanced 
against one another on the z axis, by the gravitational force and the repulsive force induced 
by the scalar charge. This balance condition is attained by choosing each scalar charge so as 
to satisfy the relationship Mj = £j(i = 1, 2), where M, is the mass of the ith singularity and 
is its scalar charge. This condition is similar to that of the MP solution, where extreme 
black holes are balanced together. In this case, the metric is described by setting k = in 
the corresponding 2-ZV solution, where Mj = ((3.10) in ||) and this spacetime becomes 
asymptotically flat. As we showed in section 2, the RLU region is not always included in 
the LU region in general, in contrast with the case with an electromagnetic field. In fact, 
the RLU region can even be seen outside the LU region (Fig.|j). However, the local matter 
effect is not strong enough to separate the RLU region from the LU region completely, so 
that chaos occurs for a particle passing through both of these regions. Hence, in this case 
we cannot clearly determine through which region chaos occurs, that is, through the LU or 
the RLU region. 



4 The shear effect of the LU region on geodesic congruence 

In section 3.2, we showed that no RLU region appears inside an LU region and K.(u,n) 
becomes negative everywhere inside the bound region in spite of the occurrence of strong 
chaos. This failure of the K,(u,n) criterion may come from the fact that the positivity of 
JC(u, n) is a sufficient condition for the local instability of geodesic but may not be necessary 
condition, because of the second term in the right side of (|2.2|) . This second term could 
induce an instability of ||n|| against the negative contribution of the JC(u,n) term. In this 
case, we could not explain chaos simply by a local instability determined from the sign of 
JC(u,n). On the other hand, the LU region still works well to judge the occurrence of chaos 
even in the Ernst case, since any chaotic bound trajectory passes through an LU region. 

In order to see the effect of an LU region even in this case, we examine the deviation of 
nearby geodesies by expansion 9 and shear of the geodesic congruence. It is well known 
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that 9 and satisfy the following equations [15 
d9 1 



dr 3 



? 2 - a, nj a^ 



-R^uW, (4.1) 



Da^ 2 



~0<V - a MP < + J V<V^ e + + -flg^ , (4-2) 



where i2^ F ^ is the trace-free part of i?^ and defined by 

4^ F) = hw h »* Rf " ~ \K^ RP ^ ( 4 - 3 ) 

where h^ v = + w^u^ is the metric of 3 space perpendicular to « M . Here we assume that 
the twist term vanishes and examine the time evolution of a geodesic congruence whose 
initial value of 9, a^ v = 0. From ( |4.1|) , if the Ricci tensor R^ v satisfies the condition 

R^u v > (4.4) 

for an arbitrary 4-velocity u, 9 negatively diverges, which is followed by the creation of a 



conjugate point of the geodesic. From Einstein's equation G^ v = SwT^, the condition (|4.4|) 
is attained, if and only if the energy-stress tensor satisfies the following strong energy 
condition (see p~5| ). 



T^uW > --T. (4.5) 



1, 
T 

As long as the condition ( |4.5D is satisfied, geodesic congruence will converge everywhere. On 
the other hand, the above energy condition tells us nothing about the shear of the geodesic 
congruence. It is possible for the geodesic to be stretched exponentially by the shear term 
a even in a spacetime with a matter field satisfying the strong energy condition (|4.5|) . This 
stretching property of shear seems to be an important cause of chaos, as we see in the famous 
baker's transformation |I7| . Since the Weyl tensor directly affects the time evolution of shear 
in (|4.2| ), it is expected to play an important role in chaos through the stretching of nearby 
geodesies in combination with their bending at the edge of the bound region. 

Especially for spacetimes with an electromagnetic field, since the time derivative of ex- 
pansion 9 is negative in (|4. 1|) , the volume of geodesic congruence always converges because of 
the strong energy condition ( |4.5| ). So shear may be one of the important factors that causes 
local instability for geodesies in spacetimes with an electromagnetic field. In the equation of 
shear ( |4.2| ), only the last two terms on the right side, 

at = ft + \k {tf) (4.6) 

are independent of the initial values of n and Dn/ Dr. Hence Q^, is important when examin- 
ing the effect of shear on chaos, since chaotic behavior of a trajectory is generally independent 
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of the initial choice of n around a given geodesic. Since each part of ( f£^ ) is trace-free, the 
sum of all eigenvalues of vanishes. Hence one of the eigenvalues of is always posi- 
tive and diverges as the geodesic approaches the gravitational source singularity, since this 
stretch comes from the tidal force by the gravitational source. 

For the asymptotically flat spacetime with a single gravitational source singularity, the 
largest eigenvalues of Q^ u , Ci, always becomes positive and decreases monotonically to zero 
as r — > oo (Fig.|5|). The direction of its corresponding eigenvector projected on the p-z plane 
is almost parallel to that of u (Fig||). (Note that those eigendirections are orthogonal to u in 
the 4-dimensional spacetime because of conditions I and II in Appendix B.) The rest of the 
eigenvalues always become negative and monotonically increase to zero as r — > oo. These 
properties for the eigenvalues of Q^, also hold for those of the Weyl tensor, since the tidal 
force in the direction of the gravitational source is characterized by the positive eigenvalue 
of the Weyl tensor. On the other hand, in multi-black holes spacetimes or Ernst spacetime, 
one of the negative eigenvalues of the Weyl tensor becomes positive inside the LU region. 
(In fact, the sign change of this eigenvalue defines the LU region.) Since is defined by the 
Weyl tensor ( |3.19| ), one of the eigenvalues of also becomes positive in the LU region. In 
the non- vacuum case, the local matter part ^R^ TF ' should also be considered in VL^. In the 
Ernst case, the local matter effect ^R^ TF ^ in fl4.6| ) is not strong enough to change the sign 
of the two positive eigenvalues, Ci and ( 2 of fll, in the LU region back to negative (Fig.|7p. 
Moreover, the eigendirections of are almost the same as those of and the Weyl tensor. 

In order to see the effect of its two positive eigenvalues on shear, we compare the eigendi- 
rections of with those of a^. As we can see in Fig.§, these eigendirections are almost the 
same, which means that the two positive eigenvalues of the Weyl tensor have the effect of 
stretching the geodesic congruence in two independent directions through the shear effect. 
As we can see in Fig.|| the second eigendirection whose eigenvalue changes from negative to 
positive is not parallel to u on the p-z plane, in contrast with the first eigendirection of £i- It 
is suggested that the second positive eigenvalue of the Weyl tensor contributes to stretching 
the congruence in the direction independent of u and causes chaos for bound orbits passing 
through the LU region. 

To substantiate this speculation, we examine the magnitude of the positive eigenvalue (2 
of inside a bound region and compare it with the value of the Lyapunov exponent A of 
chaotic motion inside this region. From Fig|](a), we can see that for a given value of L, the 
peak of (2 becomes higher and higher, as E gets larger. This is consistent with the result 
in our previous paper || that the Lyapunov exponent A, that is, the strength of chaos, is 
proportional to v*, since the effect on shear becomes larger as t>* increases. For concrete 
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values of A, in the Ernst case with (E 2 , L) = (29.12 /i 2 , 25.0 fiM), we find the numerical 
value of A ~ 3.15 x 10~ 2 /M, regardless of initial conditions, which is the same order of 
magnitude as the peak value of £2- 

As for the 6* dependence of £2, its peak becomes larger and larger as the geodesic crosses 
the equatorial plane orthogonally, although the change is not so large as to change the order 
of the value (Figj^(b)). We also find the same results in the relationship between shear 
and the LU region for multi-black hole systems. Thus these results substantiate our claim 
that the shear effect in the eigendirection corresponding to (2 is correlated with the chaotic 
motion of the geodesic. 



5 Concluding Remarks 

In this paper, we have examined the LU region criterion for chaos in the non-vacuum case, 
by analyzing test particle motion in several axisymmetric static spacetimes as well as in the 
vacuum case. 

In multi-black hole systems, the RLU region was always seen inside the LU region, while 
in a spacetime with 2 scalar-charged singularities, the RLU region was seen even outside the 
LU region. In both cases, we found that passing through an LU or RLU region gives us 
good information about the occurrence of chaos, although it is indeterminate which region 
is more fundamental for the chaos because of the slight separation of these regions. 

In an Ernst universe, however, the RLU region cannot be seen, because one of the two 
positive eigenvalues of the Weyl tensor is made negative by the strong local matter effect. In 
this case the local instability of geodesies is not explained by the curvature term in (|2.1[) or 
(2.2), since all of the eigenvalues of TZ^ become negative, which means that K.(u, n) becomes 



negative for any direction of (u, n) at any point in bound region. Then we examined the 
character of the LU region to explain such a chaotic behavior by dividing the matrix TZ^ V 
into the shear part, and the expansion part R^u^u". We find that the eigenvalues of 
are strongly affected by those of the Weyl tensor, so that in addition to the positive 
eigenvalue determined by the tidal force of gravitation, the second eigenvalue of Vt^ u also 
becomes positive around the LU region. It follows from this property that two eigenvalues 
of the shear matrix a^ u also become positive in almost the same directions as the eigen 
directions of as we numerically showed in sec. 4. Thus the LU region is characterized as 
the region in which the geodesic congruence is stretched in the direction independent of u in 
addition to the direction of u. We can speculate that this additional stretch helps to cause 
strong chaos for geodesies passing through the LU region, especially around UPO. In fact we 
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found that the Lyapunov exponent of those chaotic motions can be estimated by the peak 
of the additional second type eigenvalue. These results suggest that the Weyl tensor plays a 
much more important role in determining the chaos of geodesies than the sign of sectional 
curvature JC(u,n). 

Strictly speaking, the chaos for the orbit passing through the LU region cannot completely 
be determined by the shear effect. For example, we cannot explain, just by the LU region 
criterion, why the chaos occurs only for the orbit whose energy is almost near to Ejjpq and 
why non-chaotic orbit coexists with chaotic orbit even if both of them pass through LU 
region for the fixed value of E and L (see, for example, Fig.2 in As Vieira and Letelier 



showed in [[I4J], the chaos also seems correlated with the homoclinic tangle around the UPO. 
It might be because the stretching property of the LU region helps homoclinic mixing and 
makes the occurrence of chaos much easier around UPO. Hence, by our criterion we can just 
give a good tool to find a chaotic orbit. 

Another interesting point we wish to stress is as follows. As we showed in sec. 3.1, such a 
chaotic behavior can also be seen for the orbit with the energy beyond Ejjpq-, which cannot 
be explained by the homoclinic tangle. In a strict sense, we cannot define the chaos in such 
an unbounded case, since the orbits remain within the bound region only in the finite time 
interval. However, this result suggests that the LU region or the eigenvalues of the Weyl 
tensor could have some correlation with the unbounded chaotic type phenomena such as 
chaotic scattering or a fractal basin | T2||p2 |. Further analysis alomg this line will be left for 
the future work. 

Here we restricted ourselves to the static case, where there is no rotation effect of space- 
time involved in the geodesic equation. However, it will be interesting to extend our analysis 
to rotating spacetimes not only for the free particle motion but also for spinning particles 
2q| , because the spin-orbit or spin-spin interaction may provide new physical ingredient. 



All of these speculations lead us to future works which may reveal some important roles for 
chaos in realistic relativistic astronomical phenomena in our universe. 
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A Appendix A 



Starting from the Einstein's equation G pu = 8iiT piJ , we rewrite Eq.([2.3|) using the energy- 



stress tensor T^ u as follows 



Rfiupa = C^ vp(J + 47r{g M [ p T (7 ] ; , - g^pT^} + —Tg^ p g^ v . (Al) 



Then, we obtain 



8ir 

R(u, n, u, n) = C(u, n, u, n) + T(n, n) + ||n|| 2 T(u, u)} H T||n|| 2 . (A2) 

3 



With the definition fC(u,n) = — i?(tt, n, w, n)/||n|| 2 and fC c (u,n) = — n, u, n)/||n| 
Eq.(gg[) is now 



n) = £ c (u, n) + 4vr <| - T(u, u) \ - ^ T, (A3) 



where T(u, u) represents the energy density po as measured by the observer whose 4- velocity 
is u and T(n,n)/\\n\\ 2 represents the principal pressure p n in the direction of Then 
( |A.3j ) becomes 

K.(u, n) = K c (u, n) - 4tt{po - Pn} - — T (A4) 

Note that it is not guaranteed that the space-space components of T pv are always diagonalized 
[?]. However, it is easily shown that as long as the Riemann tensor is diagonalized, its 
components are inevitably diagonalized by the principal pressure pt for the tetrad basis 
{ e (o)> e (i); e (2), e (3)}- Hence by contracting two of the tetrad basis elements {e( ), e^, e( 2 ), e( 3 )} 
with Eq. (|Al]) , we find the relation between the eigenvalues of the Weyl tensor and those of 
Riemann tensor as follows: 

K oi = «w - Mpo - Pi} - y t 

k% = K c jk + 47i{p - Pi }-^-T (i,j,k = l~3), (A5) 

where po = ^~( e (o)> e(o)) andp^ = T(e^), &{{))■ By substituting the connection T = — po+J2iPi, 
these components become (|2.12 ). 



B Appendix B 



We first pay attention to the following properties of matrices, 7Z and C in (|3.19|) . 



I ) One of the eigenvalue of matrix TZ ( or C) is trivial, i.e., the corresponding eigenvector 
is the 4-velocity u and its eigenvalue is zero. 
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II ) Any two eigenvectors of matrix TZ ( or C) are orthogonal each other. This is not so 
trivial, since matrix the TZ is not symmetric on all indices. 

II can be proven by taking the property of the Riemann tensor and the definition of TZ into 
account. 

From Property II, it is always possible to define an orthonormal triad basis, E? (i = 
1, . . . ,3), orthogonal to u by using the normalized eigenvectors of TZ and expand by them 
the normalized deviation vector n defined as n = n/\\n\\ like 



n = 

i=l 



£n% (R6) 



from the condition, (u, n) = 0. 

From ( |£>.6| ) and the definition (|2^ ), JC(u,n) can be described as 



IC(u,n) = —R(u,n,u,n) 

i=l j=l 

3 3 
i=l j=l 



From the condition TZE^ = ot^j, ( |ij.7| ) becomes 



3 3 
i=\ j=l 



From the condition (E$, Ej) = 5ij, we can get ( |3.20|) . 
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Figure 1: (a) LU and RLU region (lightly shaded) of the 2-black hole spacetime with two 
equal masses, M, located at ±2GM/c 2 on the z axis (black dots) and the bound region 
(dotted) of a test particle with the angular momentum L = 3.4 G/j,M/c and energy, E 2 = 
0.664 (fic 2 ) 2 corresponding to -Eupo- (b)Poincare map for the orbit in the bound region with 
the initial condition p p = 0.0, z = and p = 3.0 GM/c 2 



Figure 2: (a),(b) classically admitted region in 3 black hole case with the mass ratio of 
central black hole to each of outside black hole (a) 1:0.9, (b) 1:0.95. The parameter for each 
region is E 2 = (a)0.5712, (6)0.72 and L = (a)5.2, (6)6.0, respectively, (c) Poincare map for 
the orbit in the bound region of (a) with the initial condition of each orbit is p p = 0.0, z = 
and p = 1-5, 1.6, 1.7, 1.8, 1.9 GM/c 2 . Chaos cannot be seen in this case, (d) the orbit in the 
admitted region of (b) with the initial condition, p p = 0.0, £o = and p = 1.5 GM/c 2 and 
(e) it's Poincare map. The Poincare map of the orbit passing through an LU region certainly 
becomes chaotic, before it falls into the central black hole through tiny throat. 



Figure 3: LU and bound region in Ernst case with mag netic field B = 0.15 c 2 /GM around 
Schwarzschild black hole. The parameter for bound region is L = 25.0 GpM/c and energy, 
E 2 = 29.116816 (pc 2 ) 2 . 
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Figure 4: (a) LU and RLU region in the spacetime with scalar field, where two singularities 
with scalar charge Q is fixed on z axis at z ± 2.0. The bound region (dotted) with L = 
6.9 GpM/c and E 2 = 0.90913 (pc 2 ) 2 corresponding to £"upo is also depicted. RLU region is 
almost overlapped with LU region, but not included in it. (b)Poincare map for the orbit in 
the bound region with the initial condition p p = 1.0, zq = and po = 2.8 GM/c 2 . 



Figure 5: Ci an d C2 of RN solution on equatorial plane. Here we used the parameter 0* = 0, 
L = 2.92 and (E/pc 2 ) 2 = 0.8663 corresponding to .Eupo- Ci decreases monotonically as the 
coordinate p increases, while £2 increases monotonically and approaches to zero. 



Figure 6: The distribution of eigenvector field on p-z plane for the Ernst case, Fig. 3 with 
6* = 0. The eigenvector fields corresponding to (a) (j and (b) (2 are depicted 



Figure 7: (a)£i and (2 of for Ernst solution on equatorial plane. Here we used the 
parameter 6*=0, L = 25.0 GpM/c and (E/pc 2 ) 2 = 29.116816 corresponding to -Eupo- 
(b)Ci and (2 of fi^ for the same condition as (a). Both of the cases show that Ci decreases 
monotonically, while (2 has the peak around LU region. 



Figure 8: The direction of eigenvector on p-z plane for matrix <r^ and Vt^ u along the orbit 
with the initial condition, p p = 0.0, z = and p — 4.0 GM/c 2 in the case of Fig.3. (a) and 
(b) are the two eigendirections for two positive eigenvalues of and Q^, respectively, (c) 
and (d) are the two eigendirections for two negative eigenvalues of and fl^, respectively. 
The tendency of the changes of positive eigendirections of a^ v on p-z plane along the orbit is 
almost the same as those of Q^. So is the tendency of the changes of negative eigendirections. 



Figure 9: E 2 and O* dependence of (2 of matrix C on equatorial plane for the same situation 
as fig.|. (a)£ 2 dependence of C2 for 0*=O, L = 25.0 GpM/c and (E/pc 2 ) 2 =29.116816[(i)], 
27.0[(ii)],26.0[(iii)]. We also added ( 2 of matrix C in Schwarzschild spacetime for comparison, 
where we used the parameter 6*=0, L=3.55 and {E/pc 2 ) 2 =0.9025 corresponding to -Eupo- 
(b)e„ dependence of C2 for 6* = 0[(i)], 0.0 [(ii)],7r/6.0[(iii)]7r/3.0, L = 25.0 GpM/c and 
(E/pc 2 ) 2 = 29.116816. As the 6* increases and the direction become parallel to equatorial 
plane, the eigenvalue becomes smaller and smaller. 
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